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We study several parasitic effects on the implementation of a Josephson radiation comb generator (JRCG) 
based on a dc superconducting quantum interference device (SQUID) driven by an external magnetic field. This 
system can be used as a radiation generator similarly to what is done in optics and metrology, and allows one to 
generate up to several hundreds of harmonics of the driving frequency. First we take into account how assuming 
a finite loop geometrical inductance and junction capacitance in each SQUID may alter the operation of this 
device. Then, we estimate the effect of imperfections in the fabrication of an array of SQUIDs, which is an 
unavoidable source of errors in practical situations. We show that the role of the junction capacitance is in 
general negligible, whereas the geometrical inductance has a beneficial effect on the performance of the device. 
The errors on the areas and junction resistance asymmetries may deteriorate the performance, but their effect 
can be limited up to a large extent with a suitable choice of fabrication parameters. 
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I. INTRODUCTION 

Over the last decade, important advancements in the field 
of optical frequency combs have been reported^"^. These 
have led to remarkable progresses in extending the accuracy of 
atomic clocks to the optical frequency region, with profound 
implications in several research areas, spanning from optical 
metrology^ and high precision spectroscopy^’^ to telecommu¬ 
nication technologies ^ . 

In two recent papersthe implementation of radiation 
comb generators using dc superconducting quantum interfer¬ 
ence devices (SQUIDs) or extended Josephson junctions were 
discussed. Assuming realistic experimental parameters, it was 
shown that such devices would be able to generate hundreds of 
harmonics of the driving frequency. For example, at 200 GHz 
a substantial output power of the order of a fraction of nW 
could be delivered using a standard I GHz frequency drive. 
This extraordinary frequency up-conversion opens the way to 
many applications from low-temperature microwave electron¬ 
ics to on-chip sub-millimeter wave generation. The devices 
discussed in Refs.^^’^^ were “ideal” in the sense that parasitic 
effects which can be present in a real structure were neglected. 
In light of a realistic implementation, such effects are unavoid¬ 
able and must be taken into account. In the present work we 
investigate extensively the impact of several parasitic effects 
on the phenomenology and performance of the SQUID-based 
radiation comb generator theoretically proposed in Ref. 10. 
Namely, we analyze the case in which the SQUIDs have a 
finite loop geometrical inductance and junction capacitance, 
and then we estimate the role of adding uncertainty in the 
SQUIDs areas and asymmetry parameters when building up a 
chain. We treat each one of these effects separately in order to 
emphasize their impact both on the physics and on the perfor¬ 
mance of our device. In particular we show that the junction 
capacitance plays a negligible role for our choice of parame¬ 
ters, whereas the loop geometrical inductance has a beneficial 
effect on the performance of the device. On the other hand, 
the errors on the SQUID areas and junction resistance asym¬ 
metries may deteriorate the radiation comb generator perfor¬ 
mance, but their effect remains quite moderate if such errors 



FIG. I. (Color online) (a) Sketch of the single Josephson radiation 
comb generator, a SQUID subject to a time-dependent magnetic flux 
^e{t) which induces voltage pulses V(t) across the interferometer. 
The red regions denote the two Josephson tunnel junctions, 1b is the 
constant bias current, (pi is the superconducting phase across the i-th 
junction and S are the superconducting electrodes, (b) RCSJ model 
circuit where R, Lj, and C are the resistance, the Josephson induc¬ 
tance and the capacitance of the SQUID, respectively, (c) Sketch of a 
linear array of N SQUIDs, connected together via a superconducting 
wire, and coupled to a load resistance R^. Each SQUID is pierced by 
a uniform magnetic fiux 0. The total voltage Vtot(0 which develops 
across the array is given by the sum of all the voltage drops across 
each single SQUID. 


are within a tolerance of 1% and 0.5% for the areas and the 
junction resistance asymmetry parameters, respectively. 

The paper is structured as follows: First, we review the de¬ 
vice theoretical analysis in Sec. II. In Sec. Ill we discuss how 
each parasitic effect alter the device performance: The role 
of a finite SQUID geometrical inductance and junction capac- 
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itance are investigated in Secs. Ill A and IIIB, respectively. 
Then, in Sec. IIIC we estimate the impact of an uncertainty 
in the SQUIDs areas when adding them in series to build a 
linear array, whereas in Sec. HID we consider SQUIDs with 
different asymmetry parameters. A discussion about the ex¬ 
perimental feasibility of the proposed system as well as the 
estimate of its realistic performance when all the aforemen¬ 
tioned effects are taken into account at once are the content of 
Sec. IV. Finally, our conclusions are gathered in Sec. V. 


II. THEORETICAL ANALYSIS 

In this section we briefly review the physical arguments 
leading to the prediction of the ;r-jumps of the superconduct¬ 
ing phase, and the consequent generation of voltage pulses, 
using SQUID devices. Since these were extensively discussed 
in Refs. 10 and 11 both for devices based on SQUIDs and on 
extended Josephson junctions, we recall here only the basic 
principles, without focusing on the details. 


A. Single SQUID 

We consider a SQUID biased by a constant current 1b and 
driven by an external, time-dependent magnetic flux d>e(Q 
[see Fig. 1(a)]. Due to the first Josephson relation^^, the 
Josephson current through the SQUID is 

7/ =/cisin(pi +/c 2 sin(p 2 , ( 1 ) 

where (pt and Id (/= 1 , 2 ) are the phase across and the critical 
current of the i-th junction, respectively. In the limit of negli¬ 
gible inductance^ , by introducing the superconducting phase 
across the SQUID cp = {(pi(p 2 )/2 and using the flux quan¬ 
tization relation^’^^ 92 — = — 27r<I>e(Q/<I>o, the current (//) 

vs phase relation of the SQUID can be written as 

Ij{(p’,(j)) =7+ [cos 0 sin 9 + r sin 0 cos 9 ], ( 2 ) 

where (j) = 7t^e/^Q (^0 = hjle ^2x 10“^^ Wb is the flux 
quantum), 7+ = Id Id, and r = {Id -Id)/{Ici +^ 2 ) ex¬ 
presses the degree of asymmetry of the interferometer. Equa¬ 
tion (2) describes the well-known oscillations of the SQUID 
critical current 7^(0) = maX(p7/(9; 0 ) as a function of the mag¬ 
netic flux^^, with minima occurring at integer multiples of 
<^> 0 / 2 . 

For a fixed bias current, when <I>g crosses a critical-current 
minimum we see from Eq. (2) that a change of sign in cos (j) 
must be accompanied by a change of sign in sin 9 in order 
for the current to maintain its direction. This is accomplished 
by a phase jump of ;z:i0did3-i5 owing to the second 

Josephson relation^^, results in a voltage pulse V{t) across 
the SQUID. The physical origin of the ;r-jump of the super¬ 
conducting phase can be also easily understood on an ener¬ 
getic ground. Eor a symmetric (r = 0) SQUID in the absence 
of any bias current, the time-dependent Josephson potential 
is Ej(t) = / IjV{t)dt = — EJOf{t) cos (p, where f{t) = cos(0) 


and Ejo = <^lva/{2nR). At the initial time (t = 0) this po¬ 
tential has minima at 9 = 2kn (with k integer). When the 
magnetic flux reaches the diffraction node at <I>g = <I>o/2, 
vanishes, and for > ^o/2f{t) changes its sign. The former 
equilibrium points 9 = 2k7r have become unstable and hence, 
to remain in a minimum energy state, cos 9 must change sign, 
meaning that the superconducting phase must undergo a n- 
jump to reach a new minimum at 9 = (2k + l);r. Notice that a 
finite bias current Ib is then necessary to induce a preferential 
direction to the phase jumps. 

To determine the details of the voltage pulses, we rely 
on the so-called resistively and capacitively shunted Joseph¬ 
son junction (RCSJ) model^^’^^ adapted to a SQUID [see 
Eig. 1(b)], in which each Josephson junction is modelized as 
a circuit with a capacitor C, a resistor R, and a non-linear 
(Josephson) inductance Lj arranged in a parallel configura¬ 
tion. We consider an external sinusoidally-driven magnetic 
flux with frequency v and amplitude e, centered in the first 
node of the interference pattern, so that 

<I)^(f) = ^[l-ecos(27rv/)]. (3) 

As a result, the magnetic flux crosses the nodes of the interfer¬ 
ence pattern at 4 = {2k -M)/4v, with k integer. The equation 
of motion for 9 can be written as: 

fiC h 

+ + (4) 

where C is the junction capacitance, R is the total shunting 
resistance of the SQUID, Ib is the external bias current and 
f{(pd) = ^/[ 9 ^ 0 (O]/^+- This equation can be expressed in 
terms of the dimensionless time T = 2nvt. Recalling hj {2e) = 
<I>o/ 2 ;r, we obtain^^ 

+ ^-«[/(«?,T)-5] = 0 , (5) 

where c = 2nRCv, a = I^R/ (d>ov), and 5 =Ib /7+ ^ 1 is the 
dimensionless bias current. 

The ability to generate a sequence of voltage pulses sug¬ 
gests an application similar to the frequency combs used in 
opticsIn this context, the most relevant feature becomes 
the sharpness of the voltage pulse, which is related to the num¬ 
ber of harmonics generated. The sharpness is essentially de¬ 
termined by the product 7+R, which in turn depends on the 
material properties of the Josephson junctions 

B. SQUID array with load resistor 

So far the analysis has been focused on the voltage pro¬ 
duced by a single JRCG in the absence of any external load. 
A quantity of experimental relevance is the extrinsic power 
that can be transferred to a load resistance Rl. Although the 
total output power provided by a single SQUID is fairly small, 
it can be boosted by using a linear array of N SQUIDs, con¬ 
nected together via a superconducting wire [see Eig. 1(c)]. A 
similar approach is used for the realization of the metrological 
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standard for voltage based on the Josephson effect^If we 
neglect the coupling among the SQUIDs via mutual induc¬ 
tance and/or cross capacitance and inductance of the super¬ 
conducting wire (see Sec. IV), the dynamics of each SQUID 
is independent from the rest of the arrayIn this case the 
total voltage produced by the chain is given by summing up 
the voltages developed across each single SQUIDs: 

N 

Vtotil) = E (6) 

i=l 


Assuming for simplicity that the N SQUIDs are identical, this 
can be rewritten as ytot (0 = ( 0 - This is the voltage which 

develops across the load resistor Rl [see Fig. 1(c)]. As a con¬ 
sequence, a current 1^ = NV{t)/RL flows across it, having de¬ 
noted Rl the (real) impedance of the load. The bias current Is 
is thus split into two parts, one entering the SQUID array, the 
other (II) flowing through the load resistor. This is accounted 
for by replacing the resistance R “seen” by each SQUID in 
Eq. (4) with an effective resistance 


^eff 


n N\-^ rrl 

\R Rl) Rl+NR' 


(7) 


and Eq. ( 8 ) becomes then 

hc h 


This effective change in the shunt resistance modifles the dy¬ 
namics of each single SQUID. In particular, being /?eff < it 
also reduces the power P = N'^V^/Rl that can be delivered 
to the load. Since for a single SQUID oc 7?eff, using 
Eq. (7) we And that P oc for N ^ Rl/R, whereas P oc N 
for N:^Rl/R. 

We will now turn to the analysis of different parasitic effects 
on the performance of this system. Each one of these effects 
will be first treated independently, in order to better emphasize 
its impact on the physics of the device. After that, we will try 
to give a more realistic estimate of the performance by consid¬ 
ering all these effects at once. All the results that follow are for 
an array of A = 50 SQUIDs made of Nb/AlOx/Nb Josephson 
junctionsunder a IGHz magnetic flux driving frequency. 
In particular, we have set the junction shunt resistance R = 20 
Ohm, the effective resistance Rq^ 1 Ohm (estimated using 
Eq. (7) for a linear array of A = 50 SQUIDs with a load resis¬ 
tance Rl=50 Ohm), the SQUID critical current /+ = 100 flA, 
the bias current Ib = 10 “^/+ and the amplitude of the mag¬ 
netic flux oscillations e = 0.9. Moreover, the typical size of 
the Nb/AlOx/Nb Josephson junctions is at most 1 /im x 1 /im, 
and the size of each SQUID is of the order of 1 jim — 10 jim. 


III. RESULTS 

A. Finite SQUID geometrical inductance 

In this section we investigate how taking into account a 
finite loop geometrical inductance of an individual SQUID 


modifles its dynamics under the effect of a time-dependent 
magnetic held. Labeling “1” and “2” the two SQUID arms, 
and denoting i\{t) and i 2 {t) the (time-dependent) currents 
through each of them, we define the total current I trough the 
SQUID and the circulating supercurrent is as^^^: 


I = h{t) -f i2{t), 

is = [h{t) - h{t)]/2. (9) 


Eor a symmetric SQUID (r = 0), the currents ii{t) and hit) 
are related to the voltage drops across the junctions and to the 
Josephson supercurrents by 

Vi 

h{t) =/osin(pi (0 + —, 

K 

Vo 

*2(0 =7osin(j!>2(r) + —, (10) 

K 

where the time dependence of the superconducting phases is 
given by the second Josephson relation d(pk/dt = {2e/h)Vk 
(k = 1 , 2 ), Vk being the electric potential difference across the 
k-th junction. Eollowing De Waele^ we neglect the dissipa¬ 
tive contribution to the circulating supercurrent, proportional 
to (V 2 — V\)/R. This is a reasonable assumption up to driving 
frequencies V of the order of the GHz^, and means that there is 
no appreciable contribution to is originating from the Lenz’s 
law of induction. Under this approximation Eqs. (9) reduce to 


I = 2/0 sin 
is = k cos 


91 + 92 
2 

91 + 92 


cos 


sm 


91 - 92 
2 

92 - 91 


R ’ 


( 11 ) 


where V{t) = [V'i(r) + V2{t)\/2 is the voltage drop generated 
across the SQUID. In writing the flux quantization^72. 


^ ^(0 


( 12 ) 


now the total magnetic flux piercing the SQUID is <I>(Q = 
<I>g(Q ^Lgis, which differs from the external (time-dependent) 
term ^e{t) because of the geometrical inductance of the loop 
Lg. Using Eqs. (11) and (12), after some straightforward alge¬ 
bra, we can express the total current through the SQUID and 
the total magnetic flux as: 


I = 2/0 sin 9 (r) cos 



R ’ 


<I>(Q =d>e(0 -T^/osin 


7r<t>{t)\ 


cos9(Q, 


(13a) 

(13b) 


where the phase (p{t) = [(pi{t) + 92 (Q ]/2 is related to the 
voltage drop across the SQUID via V{t) = {h/2e)d(p{t)/dt. 
Equation (13b) offers the following physical interpretation: 
At any instant of time t, the finite loop inductance modifles 
the external flux d>e(Q piercing the SQUID, and the resulting 
total magnetic flux <I>(Q has to be evaluated self-consistently. 
Once this is done, the dynamics of the SQUID phase (p{t) [as 
well as the total voltage drop across the device V{t)] can be 
evaluated via Eq. (13a). 
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In writing V{t) as the mean voltage generated across the 
two junctions, we have implicitly assumed that the two 
SQUID arms have the same inductance. Accounting for dif¬ 
ferent arms inductances Li ^ L 2 would result in an addi¬ 
tional correction to the magnetic flux piercing the SQUID, 
which would become ^’Lgis — ocLLgI/2, where 
Lg = Li L 2 , 0Cl = {Li —L 2 )/Lg, while I and is are deflned 
by Eqs. (11). From this expression we see that unless the dif¬ 
ference between Li and L 2 is large (i.e., comparable to Lg) the 
term oc aiLgl/l is a minor correction to the magnetic flux, 
with respect to Lgig. 

In order to quantify the effect of the inductance, we have 
solved numerically the RCSJ-equation (8) for an array of 50 
symmetric (r = 0) SQUIDs made of Nb/AlOx/Nb junctions 
and computed the voltage pulses for different values of the 
loop geometrical inductance, compatible with typical SQUID 
dimensions. 

In Fig. 2 we show the effect of a flnite inductance on the 
shape of the voltage pulse generated by each SQUID of the 
chain. We notice that a geometrical inductance Lg of the or¬ 
der ~pH is a reasonably good assumption for a SQUID with 
radius r of the order ^ /im if we approximate Lg fj^r, {Iq be¬ 
ing the vacuum permeability. From the flgure we see that the 
principal effect of increasing Lg is that the voltage pulses are 
delayed with respect to the Lg = 0 case, and furthermore they 
are sharper and higher. This is a direct consequence of the 
change in the time-dependent magnetic flux profile. Indeed, 
starting at ^ = 0, it turns out that <I>(r) is initially reduced by 
virtue of the second term in Eq. (13b). This means that the 
condition at which the ;r-jump of the phase is met (that is, 
<I>(r) = <I)o/2) is verified at a later time than 4 (see Sec. II), 
and the same holds for the voltage pulse. In addition, the fact 
that the shape of ^{t) is altered from the original cosinusoidal 
profile induces a faster relaxation of the phase cp toward the 
energy minimum. As a consequence, the voltage peaks for 
flnite geometrical inductance are sharper and skewed with re¬ 
spect to the Lg = 0 case (leftmost curve in Fig. 2). This has a 
beneficial impact on the emitted radiation spectrum P{Q ^), as it 
is confirmed in Fig. 3, where we show the power generated by 
a chain of A = 50 nominally identical and symmetric SQUIDs 
made of Nb/AlOx/Nb junctions, driven by a 1 GHz oscillating 
magnetic held, for different values of the loop geometrical in¬ 
ductance Lg. As we can see, the device with Lg=10 pH is able 
to provide a power of about 0.1 nW at 20 GHz (correspond¬ 
ing to the 20-th harmonics of the driving frequency). Notice 
Anally that only the even harmonics of the driving frequency 
are shown in the power spectrum of the emitted radiation, the 
contribution of the odd ones being vanishingly small for sym- 
metric (r = 0) SQUIDs'*’. 


B. Finite SQUID junction capacitance 

In this section we investigate the effect of taking into 
account a flnite SQUIDs junction capacitance. In order 
to do this, we have solved the differential RCSJ equation 
[Eq. (8)] for the SQUID phase dynamics without neglecting 
the second-order (diffusive) term. Details on the numerical 


procedure are given in appendix B. 

In Fig. 4 we show how the typical voltage pulse generated 
by each SQUID of the chain is altered due to the effect of a 
flnite junction capacitance C. We notice that increasing C up 
to 1 pF has the only effect of making the voltage peak slightly 
skewed and sharper: This would be beneflcial in terms of out¬ 
put power. For larger values of the junctions capacitance, 
the second order term in Eq. (8) becomes more important 
and the system starts operating in the under-damped regime. 
This is evident for C = 2.5 pF (rightmost curve in Fig. 4), 
at which the voltage V{t) exhibits small oscillations before 
relaxing to zero, taking also negative values. However, all 
these effects would be relevant for large Josephson junctions, 
whereas in this work we focus rather on small Nb/AlOx/Nb 
junctions, typically characterized by a relatively low capaci¬ 
tance (C < 100 fF). In this case, we see from Fig. 4 that there is 
no appreciable difference with respect to the zero-capacitance 
case (the corresponding curves are essentially indistinguish¬ 
able). As a consequence, our device operates always in the 
over-damped regime^^’^^. According to these results, we do 
not expect any relevant modiflcations in the power spectrum 
of the emitted radiation with respect to the ideal (zero capaci¬ 
tance) case, and thus we decided not to show it. 

In addition, we have also performed numerical simulations 
taking into account the combined effect of both a flnite junc¬ 
tion capacitance and loop inductance, but we did not observe 
any relevant modiflcation with respect to the results discussed 
in this and the previous subsection III A. 


C. Uncertainty on the SQUIDs areas 

When fabricating an array of N SQUIDs, it is most unlikely 
to be able to make them all identical. Inevitable imprecisions 
in the lithographic processes imply that the SQUIDs will have 
slightly different areas. As a consequence, if the array is em- 



FIG. 2. (Color online) Behavior of a typical voltage pulse generated 
by each symmetric Nb/AlOx/Nb SQUID (r = 0) of the array, for 
different values of its (geometrical) inductance Lg. The driving fre¬ 
quency is v=l GHz, whereas the other parameters are those typical 
of a Nb/AlOx/Nb Josephson junction^^, given at the end of Sec. II. 
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FIG. 3. (Color online) Power spectrum of the SQUID-based radiation 
comb generator over a 50 Ohm transmission line. Different symbols 
correspond to different values of the geometrical inductance Lg. The 
calculation is performed for a = 50 chain of nominally identical 
and symmetric Nb/AlOx/Nb SQUIDs, subject to a v = 1 GHz driv¬ 
ing. The parameters are the same as in Fig. 2. Notice that only the 
even harmonics of the driving frequency are shown, the contribution 
associated to the odd ones being vanishingly small. 



FIG. 4. (Color online) Behavior of a typical voltage pulse generated 
by each symmetric Nb/AlOx/Nb SQUID (r = 0) of the array, for 
different values of the junction capacitance C. The driving frequency 
is v=l GHz, whereas the SQUIDs parameters are the same as in the 
previous figures. 


bedded in a coil which generates an ideally uniform magnetic 
field, the resulting flux piercing each SQUID of the array 
will be different: Larger SQUIDs will be pierced by a larger 
magnetic flux, and vice-versa. This will induce a shift in the 
time at which the condition /2 (when the supercon¬ 

ducting phase experiences a ;r-jump) is met: The phase will 
jump earlier in larger SQUIDs. 

To better quantify this effect, let us associate a gaussian 
statistical distribution for the SQUID areas: 

a=Ao(i + Ca) with y{(;A) = (~^) ’ 

^ (14) 

where I^a is a dimensionless parameter quantifying the de¬ 


gree of uncertainty on the SQUIDs areas, being normally dis¬ 
tributed around zero with variance cr^, whereas Aq is the refer¬ 
ence value for the surface delimited by the SQUID loop. The 
standard deviation Ga can thus be seen as the percentage error 
within which the value of the area is known. We can write the 
external magnetic flux as: 

^e(t) = B{t)A = {Bq-Bi cos(2;rvr))Ao(l + Ct) = 

= BoAo (^1 - ^cos(27rvt)^ (1 + ^^) = 

= ^(1 + Ci) [l-ecos(2;rvr)], (15) 

where we deflned = ^o/2 and £ The phase 

jump occurs at ^e{t) = <I>o/2, that is, at a switch time t deter¬ 
mined by: 

Cl - e(l + CA)cos(2;rvf) = 0, 

—V t =-circcos (-1 (Ih) 

2nv yel + Uj 

where ^ is a non-negative integer. For sufficiently small 
the above expression for t simplifies to: 

(17) 

4v^ ^ 27rv£ 2 kv£ 

where, as in Sec. II, we have defined 4 = (1 -b 2k) /4v. From 
this expression it is evident that larger SQUIDs (Ca > 0) 
switch before (J < 4), and vice-versa. Notice also that, since 
the relation between ^a and t is linear, we can understand this 
result in terms of the distribution of the switch-times ^it), 
which can be easily computed: 

with Aa = OA/2;rev. This can be interpreted by stating that 
the times t at which the phase of the SQUIDs undergo a 7t- 
jump is normally distributed around 4 with a variance Aa 
which is directly proportional to the uncertainty Oa on the 
SQUIDs areas. 

In Fig. 5 we show how a typical voltage pulse generated 
by a linear array of = 50 symmetric SQUIDs is altered by 
assuming different uncertainties Oa on the areas, up to five 
percent. By “typical” we mean that we have first computed 
Vtot(0 for a single array of SQUIDs with random areas [ac¬ 
cording to Eq. (14)], and then we have iterated this procedure 
for many realizations of the array. We have finally calculated 
the average voltage pattern, and defined it as the typical one 
(see appendix C). We notice that the main effect is that the 
voltage peaks are broadened and lowered, due to the fact that 
a certain number of SQUIDs switch before and after 4, the 
reference switch time for a SQUID of area Aq [see Eq. (17)]. 
As a consequence, the power spectrum of the emitted radia¬ 
tion is lowered, exhibiting an exponential cut-off at high fre¬ 
quency. Despite this, we notice in Eig. 6 that this reduction 
is still very moderate for an uncertainty Oa = 0.01, in which 
case the power is reduced by less than one order of magnitude 















6 


around 100 GHz (corresponding to the 100-th harmonics of 
the driving frequency), whereas it is basically unchanged at 
20 GHz. By increasing the error to ox = 0.05, on the other 
hand, the power P is reduced in a substantial way. We fi¬ 
nally note that, in contrast to Fig. 3, in the power spectra for 
Oa >0.01 the non-dominant (odd) harmonics are visible (bot¬ 
tom curves). Remarkably, they show complex structure when 
increasing ox. This is evident for ox = 0.05: In this case, for 
> 40 GHz, the power associated to odd harmonics becomes 
of the same order, if not larger, than that associated to the odd 
ones. 


D. Uncertainty on the SQUIDs asymmetry parameters 

Another possible source of non-ideality in the fabrication 
of an array of SQUIDs stems from the asymmetry between 
the two Josephson junctions composing each element of the 
array. This is quantified in terms of the asymmetry parameter 
r = (/c,i — 4-,2)/(4 ,i + 4,2), as explained in Sec. II. We no¬ 
tice that assuming a statistical symmetric distribution for the 
parameter r around 0 (corresponding to an ideally symmet¬ 
ric SQUID) would be much detrimental for the device per¬ 
formance, because SQUIDs with ^c,2 > generate opposite 
voltage pulses with respect to SQUIDs with 4,2 < 4,i’ 
small bias current^^. Thus, when summing up all the pulses 
to compute the total voltage, the contributions associated to 
r > 0 would basically compensate those associated to r < 0, 
resulting in a poor performance in terms of output power. 

To overcome this problem, we assume that the SQUIDs are 
fabricated with a small preferential asymmetry ro, for instance 
4,2 <4,1, which correspond to ro > 0. We introduce a gaus- 
sian statistical distribution for the parameter r: 

r=r„ + C with = (19) 



FIG. 5. (Color online) Behavior of a typical sequence of two 
voltage pulses generated by an array of A = 50 SQUIDs made of 
Nb/AlOx/Nb junctions with areas statistically distributed according 
to Eq. (14) for different values of the standard deviation ox. The 
driving frequency is v=l GHz, whereas the SQUIDs parameters are 
the same as in the previous figures. 



FIG. 6. (Color online) Average power spectrum of the SQUID-based 
radiation comb generator over a 50 Ohm transmission line, for differ¬ 
ent values of the standard deviation ox of the areas distributions. The 
calculation is performed for a A = 50 chain of Nb/AlOx/Nb SQUIDs 
subject to a v=l GHz driving. 


where is a dimensionless parameter which quantifies the 
uncertainty on the SQUIDs asymmetry, being normally dis¬ 
tributed around zero with variance o^, whereas ro = 0.01 is 
the chosen reference value for the SQUIDs asymmetry. We 
have solved numerically the RCSJ dynamics of the linear ar¬ 
ray of SQUIDs following the same procedure outlined in the 
previous section. In Fig. 7 we show how the typicaf voltage 
pulses generated by an array of A = 50 SQUIDs are altered by 
assuming different uncertainties Gy on the parameter r, up to a 
standard deviation of one percent. Notice that the main qual¬ 
itative difference with respect to the previous cases, in which 
symmetric SQUIDs were considered, is that here the sequence 
of voltage pulses exhibits alternating signs. This feature was 
observed and explained in Ref. 10: Its major consequence is 
that in the power spectrum the odd harmonics are predominant 
over the even ones. 

Figure 8 shows the average power spectrum of the emit¬ 
ted radiation for an array of A = 50 SQUIDs. We notice that 
increasing the uncertainty Gy on the asymmetry parameter re¬ 
duces the power, especially at high frequency. Similarly to 
what we observed in Fig. 6, the non-dominant harmonics (in 
this case the even ones) show complex structure when increas¬ 
ing Gy. For Gy = 0.01, at high frequency (Q, > 60 GHz), the 
power associated to even harmonics becomes comparable or 
even larger than that associated to the odd ones. 


IV. EXPERIMENTAL FEASIBILITY 

In this final section we discuss the experimental feasibility 
of the setup, and we estimate its realistic performance when 
all the parasitic effects studied so far are taken into account at 
once. Some of the effects we are going to discuss were studied 
in Ref. 10, so here we just review them briefly. 

First of all, in our analysis so far we have neglected the cou¬ 
pling between the SQUIDs via mutual inductance and/or cross 
capacitance and inductance of the superconducting wire. This 
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condition, which basically relies only on the current conser¬ 
vation through each SQUID in the chain^^ implies that the 
dynamics of each SQUID is independent from the rest of the 
array, and it can be realized in practice by a suitable design 
choice. As a consequence, the voltage at the extremes of the 
array scales as the number N of SQUIDs. Accordingly, the 
intrinsic power generated by the device (that is, the power 
delivered to an ideally infinite load) scales as N^. On the 
other hand, the extrinsic power depends on the detection sys¬ 
tem used. In our case the JRCG array is supposed to be at¬ 
tached to a finite load, which effectively couples the dynam¬ 
ics of the SQUIDs: For realistic devices the extrinsic power is 
then found to scale as A, rather than (see related discus¬ 
sion in Sec. II). As shown in Sec. Ill, this A scaling is not a 
limitation in the region of tens of GHz, where sizable output 
power can be generated. Conversely at higher frequency, e.g. 
sub-millimeter region, the output power drops and the device 
design must be modified to compensate for this decrease. One 
possibility is to operate with more JRCG arrays arranged in a 
parallel configuration: In this case the contribution P of each 
JRCG array would add up and the total power would be given 
by /tot = AparF, where Apar is the number of SQUID arrays in 
parallel. 

Another important issue concerns the way the emitted ra¬ 
diation propagates across the device. When discussing the 
scaling of the power with the number A of SQUIDs in the 
chain, we have implicitly assumed such radiation to propagate 
instantly across the device. Strictly speaking, this lumped- 
element model is justified if the propagation time Tp of the 
radiation through the whole array is much shorter than the 
typical voltage pulse width, i.e. the voltage transient. This 
condition strongly depends on the specific values of the pa¬ 
rameters, which in turn are set by the device fabrication, its 
design, and the materials used. All these can be optimized 



FIG. 7. (Color online) Behavior of a sequence of two typical volt¬ 
age pulses generated by an array of A = 50 SQUIDs made of 
Nb/AlOx/Nb junctions. The SQUID chain is characterized by an 
asymmetry parameter r distribution which is gaussian and centered 
around tq = 0.01 with a standard deviation Oy [see Eq. (19)]. The 
driving frequency is v=l GHz, whereas the other SQUIDs parame¬ 
ters are the same as in the previous figures. 



FIG. 8. (Color online) Power spectrum of the SQUID-Josephson ra¬ 
diation comb generator over a 50 Ohm transmission line, for different 
values of the standard deviation Gy of the asymmetry parameter dis¬ 
tributions (centered around ro = 0.01). The calculation is performed 
for a A = 50 chain of Nb/AlOx/Nb SQUIDs subject to a v=l GHz 
driving. 


with the aim of decreasing the propagation time Xp. In any 
case, this lumped model approximation does by no means 
set any sharp boundary condition on the working operation 
of the device. Even if Xp was not much shorter, but rather 
comparable with the voltage transient, the only consequence 
would be that the interference effects shall be taken into ac¬ 
count. But since the device generates all the harmonics of 
the fundamental frequency, some of them will be partially at¬ 
tenuated because of destructive interference, some others will 
be (almost) unchanged because of constructive interference. 
Therefore the output signal may be attenuated at some spe¬ 
cific frequencies, remaining unchanged at the others, but the 
device would still work and can be used if the specific output 
frequency we want to extract has enough output power. 

In Fig. 9 we show the estimated power spectrum gener¬ 
ated by a single realistic array of A = 50 SQUIDs made of 
Nb/AlOx/Nb Josephson junctions. By “realistic” we mean 
subject to the fabrication errors discussed in Secs. IIIC and 
III D, having random (normally distributed) areas and asym¬ 
metry parameters. Furthermore, we assume them to have a fi¬ 
nite loop geometrical inductance Lg 10 pH (see Sec. Ill A). 
On the other hand, we do not consider any corrections due to 
their finite junction capacitance since we showed in Sec. IIIB 
that they were completely negligible. From the figure, we no¬ 
tice that this device is still able to provide an output power of 
about 0.1 nW around 20 GHz (corresponding to the 20-th har¬ 
monics of the driving frequency, see the corresponding black 
arrow). If we compare this to the results of Fig. 3, we notice 
that the power in this frequency range is only slightly reduced, 
as a consequence of the errors on the areas and the asym¬ 
metry parameters. A larger deterioration of the performance 
- of about two orders of magnitude - is otherwise expected 
at higher frequency (around 100 GHz, see the corresponding 
black arrow). Nevertheless, the device is still able to produce 
an output power between 0.1 and 1 pW in this range, which 
can be relevant for several applications. All these considera- 
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tions enforce the message that if the SQUIDs of the array can 
be fabricated with an accuracy of the order of 1 % on the areas 
and of 0.5% on the asymmetry between the junctions, the ex¬ 
pected performance is not altered significantly with respect to 
the ideal situation for frequencies around 20 GHz. 

Finally, we stress that all our analysis has been carried out at 
zero temperature, being more focused on the fabrication para¬ 
sitic effects. The effects of thermal noise were indeed already 
addressed in Ref. 10 for a similar setup made of yttrium bar¬ 
ium copper oxide (YBCO) Josephson junctions. In that case 
it was shown that its contribution was basically negligible, the 
signal to noise ratio being of the order of 10^ at a temperature 
of 4.2 K. Hence, we do not expect a finite temperature to alter 
significantly the results presented in this paper. 


standard load resistance of 50 Ohm. This may opens interest¬ 
ing perspectives in the realm of quantum information technol¬ 
ogy- 

The device has room for optimization by modeling the geom¬ 
etry of the single junctions, the fabrication materials, the driv¬ 
ing signal and the array design. For instance, besides SQUIDs 
made of tunneling junction considered in this work, one may 
investigate devices made of weak-link superconductor-normal 
metal-superconductor SNS junctions, such as Nb/HfTi/Nb 
Josephson junctions. 

Finally, the discussed implementation would have the advan¬ 
tage to be built on-chip and integrated in low-temperature su¬ 
perconducting microwave electronics^^"^^. 


V. CONCLUSIONS 

In summary, we have discussed extensively several para¬ 
sitic effects on the working operation of the SQUID-based ra¬ 
diation comb generator originally proposed in Ref. 10. Under 
certain conditions, we found that taking into account the finite 
loop geometrical inductance of the SQUIDs has a beneficial 
impact on the device performance, whereas the fabrication er¬ 
rors (uncertainties in the SQUIDs areas and asymmetries) tend 
to decrease it. Also, in the range of parameters considered, 
we showed that a finite junction capacitance does not alter the 
results, meaning that the device operates always in the over¬ 
damped regime. 

When all these effects are taken into account at once, we have 
estimated that a realistic array of A = 50 SQUIDs made of 
Nb/AlOx/Nb junctions is able to deliver a power of ^ 0.1 nW 
around 20 GHz, and of ~ 0.1 — 1 pW around 100 GHz, to a 
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Appendix A: Voltage spectrum and power 

To test the performance of this radiation generator, we have 
calculated the power spectrum P vs frequency O. To this goal, 
first we have computed the Eourier transform of the voltage 

V{Q.)= r dte’^’Vit). (Al) 

Jo 

The power spectral density (PSD) is then PSD{QJ) = 
\/T\V(Q)\^. Einally, the power P is calculated by integrating 
the PSD around the resonances kv (where V is the monochro¬ 
matic driving frequency) and dividing for a standard load re¬ 
sistance of 50 Ohm. This is the power we would measure at 
a given resonance frequency with a bandwidth exceeding the 
linewidth of the resonance. 


Appendix B: Solution of the second-order RCSJ equation 


FIG. 9. (Color online) Power spectrum of a single realistic sample 
of a SQUID-based radiation comb generator over a 50 Ohm trans¬ 
mission line. We have assumed a geometrical inductance L^=10 
pH, and standard deviations Oa^O.OI and (7^=0.005 on the areas and 
the asymmetry parameter, respectively (the latter centered around 
ro = 0.01). The calculation is performed for a A = 50 chain of 
Nb/AlOx/Nb SQUIDs subject to a v=l GHz driving. Blue and red 
symbols represent the odd and even harmonics, respectively, whereas 
the black arrows emphasize the frequency ranges around 20 GHz and 
100 GHz. 


To study the dynamics of the SQUID phase (p in Sec. IIIB, 
we have used a downwind finite difference approach to dis¬ 
cretize the derivatives in Eq. (8), and the resulting equation 
implemented numerically is (in the dimensionless time nota¬ 
tion): 


(p(i+i) _ 2(p(0 + (p('-i) (pW - (p('-i) 




a 


dx'^ ' dx 

cos(^^'^ sin +rsin<j)^‘^cos(p^‘^ — 5 


= 0, (Bl) 
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where is the phase at time T/, = ;r<J>(T/)/<J>o is the 

reduced flux, r is the asymmetry parameter of the SQUID, 
c = InyR^ffC is the reduced junction capacitance, a = 
/+7?eff/ (d>ov) and 5 =Ib/I^ is the dimensionless bias current. 

Appendix C: Statistical approach 

In order to estimate the effects of imperfections in the 
SQUIDs fabrication, we have followed a statistical approach. 
We describe here the procedure adopted in Sec. Ill C, the one 
in Sec. HID being equivalent. 

Given a certain value of the standard deviation ox, we have 
sampled an interval of width 8 ox by introducing a number of 
bins A^bins- We have then solved the RCSJ dynamics (8) for 
A^bins values of [corresponding to A^bins values of areas At, 
according to Eq. (14)] taken as the centers of each bin. The 
computed voltage versus time Vi{t) has been stored aside. 

At this stage, to simulate the dynamics of an array, we have 
generated N=50 values of taken from a random Gaussian 
probability distribution with zero mean and standard deviation 
ox, and to each one of these we have associated the voltage 
Vi{t) corresponding to the closest value of calculated and 


stored previously. 

For an array of N SQUIDs, under the hypothesis of inde¬ 
pendent SQUID dynamics (see Sec. IV) the total voltage is 
simply the sum of all the voltages generated by each SQUID: 

N 

Vtotit) = '£Vi{t). (Cl) 

1=1 

Indeed the presence of the load, and the fact that it effectively 
couples the dynamics of the SQUIDs, has been taken into ac¬ 
count by substituting the shunt resistance R with Rgff in the 
RCSJ equation, as discussed in Secs. II. 

Finally, this procedure has been iterated for a relatively 
large (Vi-eai=10000) number of realizations of different arrays, 
and the typical voltage of an array has been deflned as: 

'"typW = 7)^ I (C2) 

A'real j 

where the index j = 1 ... A^j-eai labels the j-th realization of an 
array. We have done this, instead of simulating the dynamics 
of all the SQUIDs of each array many times, in order to reduce 
the computational burden, otherwise enormous. 
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